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I. INTRODUCTION 



Although gross features of the string connecting a quark and an antiquark are understood, 
httle is known about the detailed structure of the flux tube. Questions remain as to the 
size and shape of the energy distribution. The linearly rising potential can be understood 
if the energy per unit length of the string is independent of the separation of the quark 
and antiquark. That still allows for a variation of the peak density and the width. We 
are also interested in the fields themselves, for example whether they are governed by dual 
superconductivity, whether they appear coulomb like at small separations as expected from 
asymptotic freedom, and what excited strings look like. 

This paper addresses the issue of the energy distribution at zero temperature. The 
difficulties of getting good signal-to-noise measurements at large qq separations are well 
known. There are trade-off's in every strategy. Our choice is to represent the qq sources 
via a Wilson loop. We then attempt to extrapolate our results to large time extent of the 
loops. The alternative is to use Polyakov lines to represent the sources. But then the 'area 
law' forces one to work at finite temperature and then extrapolate to zero temperature. 

There have been some analytical predictions of fiux tube properties. Liischer, Munster, 
and Weisz have studied the bosonic string model of the effective tube . Their results show 
that the flux tube width has a logarithmic behavior with the quark separation R. Adler 
constructed a dielectric model 0, which predicts that the energy peak density in the flux 
tube ~ and the flux tube width ~ \fR. Baker, Ball and Zachariasen have obtained the 
flux tube solution for QCD from the dual formulation of the Yang-Mills theory. They give 
analytical evidence for flux tube formation in QCD [Q. 

In this paper we report detailed calculations of the spacial proflles of energy and action 
distributions of the flux tube. Sec. || pulls together some of the simulation details for 
completeness. 

There are six orientations of plaquettes corresponding to the six components of chromo- 
electric and chromomagnetic flelds. We assign the plaquette value to the fleld-squared at 
the center of the plaquette. Hence all components are deflned at different space-time points. 
Interpolation in four dimensions is required in order to add these functions to come up with 
the energy density. The most obvious method such as cubic splines in four dimensions is not 
very satisfactory and tends to give rather obvious interpolation artifacts in surface plots. We 
decided on a method used by commercial surface plot software termed 'Kriging.' It adapted 
nicely to four dimensions without biasing any particular direction in four-space. This is 
described in Sec. |IT1[ In this way we can flU in values at every point on a lattice with lattice 
constant a/2. 

Also in this section we discuss the parametrization of the flirx tube proflle on the plane 
perpendicular to the line joining the quark and antiquark at the midpoint between them. 
It is very difficult to get accurate moments of these profiles since they are dominated by 
values of the fiux at large transverse distance where data is very scant. Our compromise 
was to find a good fitting function with a very economical parametrization. We settled on 
a simple function that is exponential for large distance and gaussian for small distances. It 
did remarkably well in that the fit was very insensitive to the sample of points used to fit 
over a very wide range of transverse distances. 

Having determined the profile parameters for a wide range of Wilson loops, we then 
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extrapolate these parameters to large time extent of the Wilson loops. In Sec. |V| we carry 
out this extrapolation and present results. 



II. SIMULATION 



The lattice observable needed to measure the flux is the following 
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where W is the Wilson loop, P!^" the plaquette located at x, P = and xr is a distant 
reference point. In the classical continuum limit 

r ^-^ -l{iFn\,-vac, (2) 

where the notation (■ ■ ■)qq-vac means the difference of the average values in the qq and vacuum 
state. From now on we shall be using field components in Minkowski space notation land 
hence 

r - \{-Bl -Bl -Bl; El El E^). (3) 

Correspondence between various components and f^'^ is standard: space-space plaquettes are 
magnetic, space-time plaquettes are electric. The energy and action densities are respectively 

e=l{E^ + Bl, 

l = liE'-Bl. (4) 

Since the magnetic contribution turns out to be negative, there is a strong cancellation 
between the two terms in the energy, whereas they are enhanced in the action. One remark 
is that the measured negative magnetic energy density was obtained from the difference 
between the value of the qq state and that of the vacuum, as shown in Eq. (|^). To understand 
the physical origin of the negative B^, one need to study the dynamical properties of SU{N) 
vacuum, which has not been fully understood. 



In the earlier paper [10| we gave some of the details of the simulation and measurement 
techniques. Briefly, our lattice size is 17^ x 20. We ran three values of /3 = 2.3, 2.4 and 2.5. 
The corresponding lattics spacings a(/?) were taken to be a(2.3) = 0.171 fm, a(2.4) = 0.128 
fm, and a(2.5) = 0.089 fm. We measured Wilson loops for all sizes up to 7 x 9 which are 
reported in Ref. [jl0| as well as the eigenvalues of the transfer matrix calculated from these 
loops. 

An important ingredient in our simulation is the analytic multihit technique in which 
one evaluates link integrals of the operators being measured. When this is applied to a single 
link on a Wilson loop the result can be given in closed form. 
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where K is the sum of six 'staples' couphng to given U in the action |TT[, and is 
the modified Bessel function. The sum of SU (2) matrices in the j = \ representation is a 
multiple of an SU{2) matrix which we denote by V 

K = bV; b= (detK)^, (6) 

Effectively what has happened is that the link has been replaced by a sum of staples asso- 
ciated with that link times a gauge invariant weight given by the ratio of Modified Bessel 
functions. 

Fig. shows the links involved. The original links in the Wilson loop and plaquette are 
the bold lines. The staple links are also shown. It is clear from this picture that not all 
the links on the Wilson loop can be evaluated with Eqn. (^. The criterion to be able to 
use Eqn. (||) to evaluate a particular (bold) link is that the corresponding staples do not 
involve any of the other bold links. A straight line of links satisfies this criterion, but links 
forming a corner do not. Therefore the plaquette and the corners of the Wilson loop have to 
be treated separately. Further since we are trying to measure the self energy of the quarks 
we need to bring the plaquette close to and in fact touching the Wilson loop which will 
result in further special cases. It is expedient to measure the Wilson loop and plaquette for 
each position and orientation independently and then evaluate the cross correlation. This 
is fine as long as a bold link on the Wilson loop does not overlap with a staple link on the 
multihit plaquette or vice versa. Since we make use of fast fourier transforms to do the cross 
correlations all relative positions of the two are automatically obtained. We must then go 
back and recalculate the correct exceptions separately. We have the choice of either doing 
more integrals analycially or drop back to fewer analytic integrations for these exceptions. 
We chose the latter, maintaining all possible analytic integrations that could be constructed 
with the above three analytic results. 

Our compromise in all these special cases was the following: We calculated (a) the 
plaquette and multihit plaquette, (b) the Wilson loop with all multihit links and the Wilson 
loop with a gap as shown in Fig. 0. Then the correlations were done using three cases: 
(i) For large separations of the Wilson loop and plaquette, all links were multihit. (ii) As 
the plaquette is brought close to the Wilson loop, the multihit plaquette was dropped, (iii) 
When the plaquette was brought in contact with the Wilson loop, certain multihit links 
were also dropped as indicated by the gap in Fig. This compromise means that having 
chosen the gap in one particular direction we can not get self energy data at other places on 
the loop for example to explore the self energy for loop with T and R reversed. 

The plaquette and corner integrals can not be done in closed form but can be evalu- 



ated using a character expansion. Details can be found in the appendix of Ref. [Q. For 
completeness we would like to pull together just the results. 

A. Character Expansions 

The basic technique used here is to expand integrands using group characters as basis 
functions. The character is given by the trace of the (2j + 1) dimensional rotation matrices. 
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(i) 



<-'-?«-^^^^^^ 0,1,1,2,.... (7) 



We parametrize group elements by an axis of rotation n and an angle of rotation about that 
axis denoted here by 2%l). 

U = cos {ip) + i sin {ilj)h ■ r, (8) 

and the group manifold is then the hypersurface of the 4-sphere and the invariant group 
integration measure is uniform on this manifold. 

We will use the character expansion of a link or product of links in the action: 

e^-(^) = T.cAP)x^HU)- cM = {4j + 2)^^, (9) 

where /2j+i(/?) is the modified Bessel function. The sum is over all representations as 
indicated in Eq. (|^) . All integrals can be done making use of the orthogonality of the group 
integration, 

/ [dU]Dlii{UyDl{UU) = ^^S,,,6^m'6nn'. (lO) 

From this relation we can find for example, 

27^2 

[dU]x^^\VU)x^^'\U) = ^-^Sn'X^'\V), (11) 

where V is also a link or product of links. This formula is very useful in evaluating the 
following integrals. 

B. Corner Integral 

Consider the integral over two links that form a corner 

{U2U1) = ^ J[dUidU2]U2Ui e-^ 

Z = f [dUM] e"^, (12) 



where the relevent terms in the action are 



S = ^(triW^U^U^) + hMUM) + Ur{U2V^)). (13) 



Z can be evaluated immediately using Eq. (11 



^ = E ^|^c,(/?)c,(/360c,(/362)x(^')(P). (14) 



where 
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P = WW2V1. (15) 
Applying the character expansion to the three terms in the action we give the final result: 

{U2U,) = - E (^^c.(/5&i)c,(/562){W^c;.x(^)(P) + 

^ . . [V2V1 - wlx^h)(P)][(j + i)x(^--^)(P) - jx(^-+^)(P)] 

where c' = dcj{[3)/d[3. 



C. Plaquette Integral 

In the plaquette integral, at most one link can be evaluated with Eq. Further it is 
clear that only two of the links forming a corner can make use of Eq. (0). Since the plaquette 
is a gauge invariant trace of now four links, this case is in fact considerably simpler than the 
corner integral. Again the details can be found in Ref. [M. Consider the following four- link 
integral 

Z{-f,(3) = J[dUidU2dU3dU4]e-^, (17) 

where 

-S^ iltrpMu^U,] + I X: tr[U,Kl]). (18) 

k=l 

We have displayed explicitly the dependence of the action on the links making up a 
particular plaquette. Each matrix is the sum of five staples. The results are 

Z{l,f3) = E9(7)(|^ n c,(/5&,)x^^'^(K.V3V (19) 

Pulhng this together, we get the result 

dZ{-y,l3) I 

/ lUdUijUriUlulU^U,) e-' = |^ J [UdU,] e"^. (20) 
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Using the plaquette and corner integrals one can reduce statistical fluctuations in the 
simulations significantly, because the measured values of plaquette and corner of Wilson 
loop involve many more links. 



III. DATA INTERPOLATION AND FITTING 

The raw fiux data can be calculated from Eq. (|Tp by Monte Carlo simulations on the 
lattice. To obtain a smooth fiux distribution and extract physical parameters, we used the 
following methods to interpolate the raw data and fit the raw data with the profile function. 
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A. Interpolation 



The six components of the flux f^'^{x) are deflned at the center of the corresponding 
plaquette and hence are all at different space-time points. Finding a good interpolation 
method in four dimensions is not as simple as it might seem. In order to combine them 
we employed a four dimensional interpolation method that minimizes the variance of the 
estimated value of the function. In other words as one varies the value of each input point 
over say one standard deviation, one requires that the linear interpolated function be as 
insensitive as possible. There is much literature on this under the name of Kriging |jl5| in 
mining engineering and a two dimensional version is employed in surface plot software. 

Consider a function, /(x), that is given as a statistical variable at a set of points, Xj. 
Suppose the mean, {f{xj)), and variance a"^. = {{f{xj) — {fixj)))"^) are known. Deflne an 
estimate of the function, /(x), as the linear combination 

f{x) = J2a,fix,). (21) 
j 

where the sum is over a subset of the points, and the coefficients a^-'s are determined by the 
set of points {xj} and the point x. Now use this estimate to evaluate the function at one of 
the given points Xi and consider the variance of the difference of the function /(xj) and its 
estimate /(xj) 

^'e = mi^^) - (/(^.))) - - (/>.))))') 

= - 2 E «i ^x, ,x, + E % ^x, ,x., ■ (22) 
j jj' 

The sum excludes the point Xj. The covariance matrix is defined 

= (ifi^k) - (/(Xfc)))(/(x0 - UXxi)))). (23) 

The Kriging method determines the coefficients aj by minimizing subject to a con- 
straint, i.e. we minimize 

£(a,) = a,2(a,) + 2^(^a,-l), (24) 

j 

where is a lagrange multiplier. The constraint arises from the condition that when applied 
to all points in the set, /(xj) underestimates /(xj) as often as it overestimates it: 

{{J2a,fix,)-f{x,))) = 0. (25) 
j 

Using ((/(xj))) = / this gives the constraint J2j % = 1- Hence the coefficients are given by 
the equations 

j' 



E«/ = l- (26) 
i' 
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We want to apply this to interpolate f{x) at a general point. The a/s are determined 
from the covariance matrix at the general point. Since the a^'s are independent of the overall 
normalization of the covariance matrix and the relevant property is how the covariance 
matrix varies with distance between the two arguments, so in practical applications one 
assumes that the covariance matrix falls as a function of distance, going to zero at large 
separations. 

Our choice is: 

~exp(-|j-/|/r). (27) 

One can easily check that if one chose the set {xj} to include the point Xj, then the cor- 
responding aj = 1, and all other a^-'s equal zero. For a point near one of the Xj, the 
corresponding aj = 1 dominates. Since we have an ad hoc ansatz for the covariance matrix 
the actual errors in our data do not enter in the interpolation. In fact this is just a method 
of linear interpolation which weighs the nearby points more than the far points which is 
divorced from its statistical origins. 

The freedom in this method is in the choice of neighboring points and the value of the 
correlation length r. Our criterion was to get a smooth interpolation which did not produce 
variations of the order of a lattice spacing. We found r ~ the lattice spacing and a radius that 
included one or two dozen points in the four dimensional neighborhood of the interpolated 
point. The errors were calculated by linearly interpolating the errors from the neighboring 
points. This overestimates the errors since it assumes that the values used to interpolate 
are completely correlated. Assuming they are completely uncorrelated then interpolating 
the variance gives an error typically less than half of the linearly interpolated errors. Since 
we do not know the covariance matrix we took the worst case. 

Using this method we filled in every point at the half lattice spacing for each component 
which then allowed us to combine any components. Figs. § and |^ give surface plots of the 
interpolated data for a Wilson loop of T = 7 and R = 5 with (3 = 2.4. Shown are two 
sections, one through the two quarks and the second on a plane midway between the two 
quarks to show the flux tube profile. 

Fig. 1^ shows the electric and magnetic components —^Bp ^Ej_, and —jB'j_. From 
these plots we can see that only the electric components have prominent peaks around the 
quarks, and the two transverse components, ^E'j_ and —^B'^, are essentially the same. These 
two components cancel in calculating the energy and add in the action. They are the widest 
of the four profiles and hence the action has a wider profile than the energy. 

Fig. ^ shows the action and energy for the same parameters as above. The peaks around 
the quarks rise above the background equally in the energy and action since they come 
predominantly from the electric components. The Michael sum rules confirm this behavior 
as described in the earlier paper |]10|]. The flux tube profiles for action and energy are quite 



different. Note the large cancellation in the energy compared to the action and the resulting 
difficulty in measuring the former. Rotational symmetry is reasonably well restored in the 
action profiles but there are large lattice artifacts in the energy profiles. 
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B. Fit to Profile Function 



The basic parameters we are interested in is to find the peak value and widths of the 
flux density of the energy and action. Obtaining good data on the second moment of the 
profile in order to determine its width is problematic. The proflles fall sharply and the 
fourth moment needed to estimate error would be sampling only large distances from the 
axis where the signal to noise is small. Instead we fitted the energy and action density in 
the plane at the midpoint between q and q using the function 

/(rx) = aexp(-v/c2 + (rx/&)2). (28) 
or rewriting in terms of the width at half maximum, B 

f{r^) = y4e^exp(-^C2 + {2Cln2 + ln^2){r^/ B)^) . (29) 

The peak value, A, and the width B, were very well determined using a fit fo^' each 
of 70 cases of different loop sizes and values of p. (We used Minuit to determine the fit 
and errors.) For the third parameter we chose the decay length of the tail of this function, 
D = {2Cln2 + lii?2) / B. This parameter was less well determined. We attempted to 
include an additive constant to the function since our sum rule paper gave evidence for 
it. However this gave too much freedom to the fit and gave unstable answers. We also 
tried a gaussian plus background constant and this did not fit nearly as well as the function 
chosen. In conclusion, the data is not good enough to determine four parameters, and the 
best function we found has some ambiguity in the choice but the peak value. A, and width 
at half maximum B are very well determined. The exponential fall off of this function is 
perhaps modeling a true exponential falloff plus possibly a background which we believe is 
there due to the fact that we are measuring the flux relative to a reference point that is not 
at inflnity. 

Our complete fltting results are shown in Tables ||-[V^for action and energy distributions. 
The four numbers for each loop size are the peak value, A in GeV/ fm^, the width at half 
maximum, B in fm, the exponential decay length, D in fm, and the integration of the area 
under the curve to give the string tension in GeV/ fm. As a reminder, since we are looking 
at the middle time slice and looking at the mid-point of the flux tube, these parameters are 
measured in the x, y plane for z,t at the middle of the Wilson loop. 

Our characterization of a good flt requires explanation. We do not report the because 
the flt is used to average out the lattice artifacts in the rotational symmetry. For small loops 
the errors in the data are quite small but rotational symmetry is not equally good. As a 
consequence can be very large for flts since the fltting function is rotationally invariant. 
When the errors in the data approach the size of the lattice artifacts x^/degree of freedom 
drops to the order 1. We beheve the flt is good because of the stability of the numbers 
over changing the set of points used in the flt. We determined the parameters using 9 sets 
of points, choosing the interior of circles of radius 1.5a to 5.5a. (Recall our interpolation 
produces values at every point at half a lattice spacing.) For all our flts the variation in the 
parameters over this set of points was of the order or less that the statistical Minuit errors 
quoted for the parameters in the tables. We did not fold in this systematic error. However 
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we believe we have overestimated errors in the interpolation and we believe the effects of 
these two contributions are comparable. 

Fig. ^ illustrates a general feature of our data. The cluster of three points for each R 
and T correspond to the three quantities: 

e=\{El + Bl) + \{El + Bll 

e(T ^R) = \{El + El) - \{El + El). (30) 

If one turns the Wilson loop on its side, the || components are unchanged but the _L com- 
ponents of the electric and magnetic fields are reversed: E^ —E\. Hence there is a 
sign change in the third expression. The central points of the cluster are the || components 
only. The clustering of the points implies that the _L components of the electric and magnetic 
contributions to energy density are approximately equal hut of opposite sign and cancel. The 
width of the peak is even less sensitive to the transverse components giving essentially the 
same value for all three points. 



IV. RESULTS FOR THE GROUND STATE 

A. Extrapolation to Infinite Wilson Loop Time Extent 

In the companion paper we described the use of the trasfer matrix eigenvalues to 
extrapolate the data to infinite time. Similarly to the potential case, the finite time extent 
of the Wilson loop W{R,T) introduces contamination of the f^^{x) by the excited states of 
the colour field [jl6|. From the transfer matrix representation of the correlation between the 
Wilson loop W and a plaquette P, 

{WP) = Z"^Tr{r'''-^'''Sr^J^''PT^J^''S), (31) 

where 5" is the operator which excites the qq states from the vacuum, and Tqq is the transfer 
matrix in the qq sector. We obtain for any component f^" (= F) 

F = J^ + 7^le-^^^ + + ■ ■ • (32) 

where denotes true ground state average, and Ei stands for the energy gap of Tqq, which 



can be determined by fitting the Wilson loops to the exponentials as described in Ref. [10 
Formula (|32D applies to any component (/iz/) and any three-space location x of the plaquette 
P^^{x, Xq). To minimize the effect of the cross terms, Xq was always chosen at the middle of 
the time interval of the Wilson loop (hence the term e~^^"^/^). For any parameter extracted 
from the flux distribution f^'^, e.g., the peak density, one can use Eq. (|3^ ) to fit the parameter 
values corresponding to various Wilson loop time extents, and extract the asymptotic value 
T ioi T ^ oo. 

A sample of the results of extrapolating are shown in Figs. |^, ^ and 0. One can see 
that the extrapolation is reasonable. However, there are some cases in which the fit did 
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not work, one could trace the problem to too large a variation for small T where the points 
are most accurate. This indicates that a constant and one exponential is not sufficient to 
take for small T. We then dropped points as needed at the small T. These cases were all for 
f3 = 2.5 which correspond to the smallest physical times. 



B. Discussion of Scaling Plots for Parameters 

After we extrapolate the measured physical quantities to infinite Wilson loop time extent, 
we can study the behaviors of these quantities with the quark separation R. In Figs. || and ^ 
we just the results for the peak densities and flux widthes. 

The basic issue is whether the peak value of the energy density stabilizes to a constant 
or goes to zero with quark separation. This is not easy to settle as can be seen in Fig. |^ 
(a). Roughly speaking we know from the linearly rising potential that the string tension ~ 
(width)^ X (peak value) should be constant. Both Figs. ^ and ^ show that we are marginally 
asymptotic in quark separation, R. The two curves in Fig. | are ~ and ~ l/i?. A 

coulomb field would fall like the former but since the string tension is constant the asymptotic 
width would have to grow like R^ which clearly it does not. Therefore we can rule out a 
coulomb field as expected. Interestingly for small separations, the eyeball fit to ~ l/i?^ is 
quite good which may be due to a coulomb like behavior at small distances. (The above 
argument that the width must grow like R^ does not apply because there is no string for 
small R.) The dielectric model predicts the peak density ~ 1/_R. Such a behavior would 
imply the width ~ a/R which is certainly possible in our data. We can say that the peak 
energy density and width are consistent with a constant value for large quark separation 
but we can not rule out a slow variation. The issue can be tightened by making use of the 
Michael sum rules as we mention in the next section. 



C. Using Sum Rules to Tighten Predictions of Energy Behavior 

A consistency check on the flux distributions can be obtained by using the Michael sum 
rules |T^ for energy and action. 



X 



1 J:{E{xY - B{xf) = -Aeo{R) - ^] - (33) 



Here Eq{R) is given by 



(W^(i?,r)) = ^A,e-^^(^)^, (34) 



and (■ = ^). In Ref. [|T^ we have shown that our data are essentially consistent with these 
sum rules. The one difficulty is the fact that the self energy, c(/5)/a(/5), determined from 
the potential differs from the self energy determined from the the action sum rule. This 
may be due to an ambiguity in the definition of self energy or possibly due to our classical 
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expressions for energy and action which ignores quantum corrections. By taking a derivative 
of these expressions with respect to the quark separation, R, this difficulty is avoided. This 
gives the relation 

a A = -f-o- aA^\ Y.^E{xf - B(xf)- o ^\ Y.{E{xf + B{xf). (35) 

The sums are now over the plane midway between the qq pair. 

Fig. [TO shows the data for the ratio, oaI^^ which were measured for /3 = 2.3, 2.4 and 



2.5, at various qq separations R. One can see that for small R the ratio is small, however, 
at large R most data show large values (~ 10.) for this ratio, but the fluctuations at large 



R make an estimate difficult. From the sum rules above and Fig. [T^ we estimate that the 
/? function, —(3a/ a) ^ 10. ± 2., compared to the current estimates 7.7 ± 1.0. The 

asymptotic value is —51/121 + 3tt'^(3/11 ~ 6.0 for (3 = 2.4. There is ample evidence from 
other measurements that although scaling works well, asymptotic scaling is violated |l^ and 
hence we do not expect to get the asymptotic value. 

An alternative approach is to assume the sum rules are correct and use them to infer 
information about the energy density from the action density which is far easier to measure 
since relative errors are down by an order of magnitude. As is clear from the sum rule the 
action does not scale yet the variation over these values of (3 is very small. An examination 
of Fig. ^ (a) shows that the action for each (3 seems to stabilize to a constant for increasing 
distance for the peak density and for the width. This is strongly supported by data for 
P = 2.3 and 2.4. For /3 = 2.5 the quark separation R appears to be too small to draw a 
conclusion. From these data one can see that the peak values seem to approach a nonzero 
finite constant at large R. We would hke to use the sum rules to predict the behavior of 
the energy density. A constant peak energy density follows only if the widths of the energy 
and action peaks have the same behavior. Figure |TT| shows that in fact they do. From 
this and using the sum rules we conclude that the energy density stabilizes to a constant 
value also. This conclusion is an argument against the dielectric model 0. However due to 
the limitation on the available range of R of our data we can not make conclusions about 
logarithmic behavior of the flux tube width as predicted by Liischer B. 



D. Longitudinal Profile of Flux Tube 



To see the proflle of the whole flux tube, in Fig. we plot the integrations of flux action 
on various transverse slices along the qq axis, where Z = is the location of one quark 
source, the positions of Z < corresponds to slices outside the qq pair, and the positions 
with Z > corresponds to slices between the qq pair. Because of the symmetry of flux 
distribution about the midway point between the qq pair, we only show the distribution 
from one end to the midpoint between the qq pair. One can see that the peak values are in 
the neighborhood of the source, which is caused by the large self-energy concentrated around 
the source. However, the flux action values between the qq pair are comparable to the peak 
values. This implies that the integrations of flux action on the transverse slice near the qq 
pair are not dominated by the self-energy contributions, although the density values due to 
self-energy are large, their distribution regions are very narrow. This agrees with the action 



12 



density distribution shown in Fig. |^ (a), where the peaks of action density are restricted in 
a narrow region at each source. From Fig. |l^ one can also see that away from the outside 
of the qq sources the flux action decreases rapidly with the distance from the source, and 
almost vanishes at the moderate distance, R = Aa (i.e., Z = —4). However, the flux action 
between the qq pair is large, and almost keeps to be a constant as Z increases. This agrees 
with the expectation for the flux tube formation. 
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TABLES 



TABLE L The fitting parameters for the action distributions with (3 — 2.3 and various 

Wilson loop sizes, T and R. For each loop size we list in order the values of the three 
parameters, A in Gev/fm^, B in fm and D in fm. The last value in each case gives the 
integration of the area under the curve in Gev//m. 

TABLE II. The fitting parameters for the energy distributions with /3 = 2.3 and various 
Wilson loop sizes, T and R. Which is similar to TABLE I. 

TABLE III. The fitting parameters for the action distributions with /3 = 2.4 and various 
Wilson loop sizes, T and R. Which is similar to TABLE I. 

TABLE IV. The fitting parameters for the energy distributions with /3 = 2.4 and various 
Wilson loop sizes, T and R. Which is similar to TABLE I. 

TABLE V. The fitting parameters for the action distributions with j3 = 2.5 and various 
Wilson loop sizes, T and R. Which is similar to TABLE I. 

TABLE VI. The fitting parameters for the energy distributions with j3 = 2.5 and various 
Wilson loop sizes, T and R. Which is similar to TABLE I. 
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FIGURES 



FIG. 1. Links involved in link, corner and plaquette integrals. 

FIG. 2. Surface plots of electric and magnetic components. For each component two sections 
are shown, one is on the longitudinal plane containing the qq pair, another is on the transverse 
plane midway between the qq pair, (a) and (b) are for ^-Ey , (c) and (d) for — ^-By , (e) and (f) for 
^El, (g) and (h) for -^5^. The data were measured on the Wilson loop with T = 7, R = 5 and 
p = 2.4. 

FIG. 3. Surface plots of energy and action, similar to Fig. ^, and obtained from the same data 
set. (a) and (b) are for the action density, (c) and (d) for the energy density. 

FIG. 4. Peak density for R x T Wilson loop sizes with /3 = 2.4 and i? = 3 - 6, T = 3 - 6. 
For fixed T the points decrease monotonically with R. Triangles: energy density for R x T loop; 
circles: energy density for T x R loop; squares: || components of E'^ and only. 

FIG. 5. Extrapolation of action peak density to infinite Wilson loop time extent. The extrap- 
olated values are plotted at T=12. The data were measured on the Wilson loop of the size RxT 
with (3 = 2.4, for T = 3 — 9 and R = 2> (open triangles), R = A (solid triangles), R = 5 (solid 
squares), R = Q (open squares) and R = 7 (solid circles). 

FIG. 6. Extrapolation of action flux tube width at half maximum to infinite Wilson loop time 
extent, labels are similar to Fig. ^. 

FIG. 7. Extrapolation of exponential tail decay length of action flux tube to infinite Wilson 
loop time extent, labels are similar to Fig. |5|. 

FIG. 8. Peak values of energy and action density; solid squares: /3 = 2.5; triangles: f3 = 2.4; 
open squares: P = 2.3. The two curves are 1/R and arbitrarily normalized. 

FIG. 9. Width at half maximum for energy and action density, labels are similar to Fig. |8[ 

FIG. 10. Ratio of integrations of action to energy densities on center slice, labels are similar 
to Fig. |. 

FIG. 11. Ratio of widths of action to energy profiles on center slice, labels are similar to Fig. ^. 
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FIG. 12. Integrations of action density on transverse slices along the qq axis. The data were 
measured for (3 = 2.4 on the Wilson loops extrapolated for large time for, i? = 3 (open triangles), 
i? = 4 (solid triangles) , R = 5 (solid squares) and R = 6 (open squares) . 
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